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Abstract —We reduce the problem of optimal beamforming 
for two-way relay (TWR) systems with perfect channel state 
infomation (CSI) that use analog network coding (ANC) to a 
pair of algebraic equations in two variables that can be solved 
inexpensively using numerical methods. The solution has greatly 
reduced complexity compared to previous exact solutions via 
semidefinite programming (SDP). Together with the linearized 
robust solution described in [1], it provides a high-performance, 
low-complexity robust beamforming solution for 2-way relays. 

Index Terms —Two-way relay, Beamforming, Low complexity. 
Conjugate gradient. 

1. Introduction 

The analog network coding (ANC) technique has been 
proven to lead to significantly higher throughput in wireless 
router scenarios [2]. In reference [3], ANC is evaluated in the 
context of two-way relay (TWR) systems with two single¬ 
antenna source nodes communicating via a multi-antenna 
relay. In that paper, an optimal beamforming solution was 
derived that made use of the S-procedure to reduce the 
beamforming problem to a system of linear matrix inequalities, 
which can then be solved by semidefinite programming [4]. 
The paper also noted that the optimal solution could be 
expressed in terms of four complex design parameters. The 
current paper provides significant simplifications over that 
result. We reduce the number of design parameters from four 
complex parameters to two real parameters. Furthermore, we 
reduce the problem to an unconstrained minimization problem 
in two real variables, giving rise to a system of two algebraic 
equations in two unknowns that can be solved inexpensively to 
arbitrary accuracy using conjugate gradient or other numerical 
methods. 

The optimal solution described above applies to the case 
where the channel state information (CSI) is known exactly, 
which is commonly designated as the “nonrobust” case. In the 
“robust” case, the CSI is only known to a certain tolerance. 
It was shown in [1] that given an exact solution for the 
nonrobust case, a low-complexity suboptimal robust solution 
with very high performance can be found. Thus our algebraic 
nonrobust solution can be used as part of a complete low- 
complexity solution to the robust beamforming problem for 
two-way relays with ANC. 

The rest of the paper is organized as follows. In Section 
II, we present the system model and formulate the problem; 
in Section III, we reduce the problem to a much simpler 
problem; in Section IV we give algebraic solutions to the 
simpler problem; in Section V we present simulation results 
for the nonrobust case; in Section VI we describe the subop¬ 
timal robust solution, and present simulation results; and in 
Section VII we summarize our conclusions. 


The notations used in this paper are listed as follows. We 
define (-j^, {■)^, (•) as the transpose, Hermitian transpose, 
and conjugate operations, respectively. Re(-) is the real part 
and Im(-) is the imaginary part of a complex variable. We use 
Tr[-] to denote the trace of a matrix. 

II. System Model, and Statement oe the 
Beameorming Optimization Problem 

We consider a two-way relay system similar to the one 
introduced in [3], which consists of the relay node R and 
two terminal nodes S'! and 52. The relay is equipped with 
M antennas and the terminal nodes are each equipped with a 
single antenna. For terminal node Si {i = 1, 2), we define pi 
as the transmit power level and hi G as the complex 

channel gain from node to relay. We further define erf as the 
noise variance in the received signal at 5^ (i = 1, 2), and 
as the noise covariance for the received signal at R, where 
all noises are assumed to be circularly symmetric complex 
Gaussian (CSCG). It was shown in [3] that for an ANC system 
in which the terminal nodes exchange information in two 
consecutive time slots under conditions of channel reciprocity 
(justified in [5]), transmit power at the relay R is given by 

G(A) = ||A/ri||V + l|A/r2fp2+Tr[A^A]a|, (1) 

where A G is the relay’s beamforming matrix. 

Reference [3] also shows that if the SINR at node Si is 
constrained to be at least 'ji{i = 1, 2), then assuming perfect 
knowledge of CSI (which is denoted as the “nonrobust” case) 
the optimization problem to minimize the relay power can be 
formulated as follows: find (i=l, 2) 

= argmin[G(A)] s.t. /,(A) > 7,crf, (2) 
A 

where 

/i(A) = |h,f - lll^rA||V|7*,, {k = 3-i). (3) 

We note that the problem in (2) is not convex in general, 
because the constraints are not convex functions. 

III. Reduction to real-valued rank 2 problem 

In this section we show how (2) can be transformed into a 
much simpler problem with real coefficients. 

It has been shown previously in [3] that the solution A* of 
(2) is of rank 2. Specifically, A* can be expressed as 

2 

A* = {at.)ijhihf = [hi, /r2]a»[^f ; h^], (4) 

*j=i 

where a* is a complex 2x2 matrix. The objective function 
condition and constraints in (2) can be rewritten in terms of 
the matrix a*. The coefficients which appear in this simplified 



version of (2) will be complex in general; but it is possible 
to further simplify the expressions so that all coefficients are 
real as follows. First we define 

_ ||(/ri/||/ri||±e^%/|l/r2||)|| 

~ \h^h,\’ ^/2 

We then choose the following orthonormal basis {e+, e_} for 
the space spanned by hi and e^^h 2 '. 

e± = {hi/\\hi\\ ± e'’^/r2/|l/j'2||)/(v^f±). 

The following may be verified, where r = t-/t+ (note r >0): 

= 0 ; 

hi = ||h.i||f+(e+ + re_); h.2 = e'’^|jh.2||f+(e+ - re_). 

Since the vectors {e+,e_} defined above span {hi,h 2 }, we 
may alternatively write 

A* = [e+ e_]a[e^ ; 

where a is a 2 x 2 complex matrix. Using the following 
rescaled constants (i = 1, 2; fc = 3 — t) 

q, =p^\\h,ftl/a%; Ci = pfc|j/rif |jh.2f 

di = WhiW^tlayaf; n = [1; r];, T2 = [1; -r], 
the optimization problem becomes: 


a* = argmin[p(a)] s.t. /i(a) > 1, (z = 1,2), (5) 

a 

where (z = 1, 2; A: = 3 — z) 

g{a) = gillarif + q 2 \\aT 2 \\‘^ + Tr[a^ a]; 

Mo-) = - diWM “IP- 

Note that although a\g{a) gives the actual power, for conve¬ 
nience’s sake we will refer to g{a) as the “power function”. 

The functions in (6) can be compactly expressed as 
quadratic forms. First we define (z = 1, 2; fc = 3 — z) 

Tii = TiM; -m = qiTii + q 2 T 22 + I. (7) 

Next, for any 2x2 matrix a we define the operations: 


[oil 012 ' 

021 0:22]^; 


a 0 


Oiil O21I 

0 a 

; a = 

012I 022! 


Finally we define 

7VT = ZTZ; ‘^ki = Qi = d^Ta^ (9) 

where M, T^.^, and are all real symmetric 4x4 matrices. 
Using this notation, we have 

g{a) = Ma; fiM) = ol^Q^ol. (10) 

Note that all the coefficients in (10) are real. We now show that 
for any locally-optimal complex feasible solution to (5) with 
g and as in (10), there also exists a real feasible solution 
that achieves the same power. This implies there always exists 
a globally optimal real feasible solution. 

Let us write x = Re [a] and y = Im[a]. Then we may 


consider g and fi as functions of x,y: 

g{x,y) =x'^Mx + y'^My, 

Mx,y) = x'^Q^x+ y'^Q^y, (i = l,2). 

In terms of M and Q, the KKT conditions corresponding to 
the minimization problem (5) are (Ai,A2 > 0) 

Mx = \iQiX + \ 2 Q 2 X] My = \iQiy + \ 2 Q 2 y, (12) 

where Ai(l — fi{x,y)) = 0, (z = 1,2). Consider first 
solutions of (12) where fi{x,y) = 1, z = 1,2. Using (12) 
we may verify that the real beamforming matrix Ca,® + (yy 
also satisfies the KKT conditions, for any real and (y. 
Furthermore, if fi{CxX + Cyy,0) = 1 (* = 1,2) are satisfied 
with equality, then it can be shown that (^x -b Qy is a real 
solution which satisfies both constraints and has the same 
power as the complex optimal solution. The proof that suitable 
CxXy can always be found is accomplished on a case-by- 
case basis using various geometrical arguments. In the case 
where only one of the two constraints is satisfied with equality, 
similar arguments may be used. Details may be found in [6]. 


IV. Algebraic solution of the real problem 

^(1) ^(2) 

Define real 4x4 orthogonal matrices U ,U as (p = 
1 + r^) 
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(13) 




The columns of ' form an orthonormal basis of for 
A: = 1,2, so for any vector x G we may write x = 

^x^^^^ It is possible to show that {k = 1, 2) 

fk{x) = X^Q^X = p {{ckP - dk){xP)‘^ - dkixP)'^'^ . 

Consider first the case where both constraints fk(x) >1, {k = 
1,2) are satisfied with equality. Then we have: 


where 


and 


Ofc = 


x[ = ±Pk{x2 )> 

Pk{x) = \fakMbkM (fc = l,2), 

dk 


(14) 


^kP dk 


bk = 


^kP dk 


(fc = l,2). 


The other components of x^^'> are uniquely determined 
via (15), where u 


(k) 


M) 


denotes the j’th column of U 

;(2) _-(2 )t 
•^3 ) “4 ) “3 5 “4 J 

early independent when r 7b 0,00.) Using the first two 
rows of this matrix, we may construct a matrix P such 


(Note that the vectors lin- 


that xM = p 


(1) ^(1) ^(2) ^(2) 


. Defining Z = 



(15) 
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{U^ '^)'^MU^ P, we may then write the power function 
as: 


_ r^(i) ^(1) J2) ^(2)i;7r^(i) ^(i) ^(2) ^(2 )iT 

yV'^t — [•''1 J •''2 ) •''1 ) •''2 ) >''2 5 •''1 > •''2 J ■ 


The solutions we are seeking are the unconstrained solutions 
to 

Hi{x) 

X 

±Ai2(y) 

y 

(16) 

Note we have left off one of the ±’s because of symmetry- 
these other solutions will be the negatives of the solutions 
to (16). For each choice of sign in (16), it can be shown 
there is a unique solution. This can be seen geometrically 
as follows. The solution x considered as a point in R"* is 
determined by the intersection of the nested family of strictly 
convex sets {g{x) < with a strictly convex set S 

that is the intersection of two strictly convex components of 
the constraint sets (one component from each constraint in 
(5)). The set {g{x) < K} is equal to the point 0 € 
and since 0 ^ S' it follows that the minimum of g{x) on S 
is positive. Because of the convexity of the sets involved, the 
smallest value of K for which the intersection {g{x) < FfjnS 
is nonempty produces an intersection consisting of a single 
point, which is the unique global minimum of the function 
g{x) under the given constraints. 


„(1) ..(2)^ 


) = argmm 


gi[x) 

X 

±^2(2/) 
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In terms of the solution to (16), the optimized beamforming 
matrix (in vector form) is 


OL = X = U^^^P 


±g2{x^2^) 


Using the symbolic algebra software Maxima, we may find 
explicitly 



1 r 1 r ' 
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There are two local optima, corresponding to the ± sign in the 
expression. Maxima can also be used to find the expression 
(17) for Z, which utilizes the definitions 


p = + 1; p' = — 1. 


We consider now the possibility of optimal solutions for which 
one constraint is satisfied with equality and the other with strict 


inequality. First we suppose that fi{x) = 1 and f 2 {x) > 1 
for the optimal solution x. Since x satisfies fi{x) = 1, we 

may write +^1(2^44^) where x = x^^\ We 

may also write the power as g{x) = where 

Y = iU^ ^Y'MU^ 4 Since M is positive definite, then so 
is Y. It follows that for x to minimize the power subject to 
the constraint fi{x) = 1, we must have = x^/^ = 0. We 
thus have 


-1 

T 

2/11 
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Ml ( 4^4 
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where using Maxima we find 

yii = p~^ + q2P^ + p) ; 2/12 = 2/21 = p"^( 2 girp'); 

2/22 = + p). 


„(i) 


satisfies 


+ 2/22a;4^ = 0) 


This expression is minimized when X 2 

(1) , ( (l)^ , yi2bl{x^2^f 

Viibix\ + yi2Pi[x)2 ’) H-771)— 

Pl{x\ ) 

which leads to 

0 = ((2/ii('i + 2/22 )^ - 4j/i2^i) 

+ {{yiibi + 2/22)^01 - 4^1201^1) (a; 44 ^ - 2 / 120 ?- 


Notice there will always be one positive and one (nonphysical) 
negative solution for (a;44^- We should consider both the 
positive and negative solution for a;44 and we may choose 
the positive root to obtain = Pi{x^-p)- Finally, we must 
check whether or not the second constraint is satisfied: 


1 <f2{x) 

=p [{C 2 P - d 2 ){x^i^U^^^ ■ 
— dk{x^i'’U^i'^ ■ 


V(2) 


+ X 


(l)-(l) jv( 2 ) 


M, 


UX 


nr 


which is equivalent to 

^ < ((C2P- d2){p'x[^'' - 2rx44^ - d2{2rx^^^ +p'a:44^) • 

For solutions which satisfy the second constraint and not the 
first, the same equations are used except the indices 1 and 2 
are exchanged. 

In summary, we have found six candidate solutions to be 
evaluated: two which satisfy both constraints with equality, 
and four others that satisfy one constraint with equality, and 
one with strict inequality. The candidate that has the lowest 
power will be the true optimal solution. 

In practice, the candidate solutions must be computed 
numerically. Four of the candidates (those for which one of 
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the constraints is not satisfied strictly) are obtained using the 
quadratic formula. The other two (specified by (16)) can be 
evaluated to the desired precision using numerical methods 
such as conjugate gradient. Convergence of the Polak-Ribiere 
and Conjugate Descent variants of the conjugate gradient 
method is guaranteed due to the boundedness of level sets 
and Lipschitz continuity of the gradient of the function to be 
minimized. [7]. 


V. Simulation of nonrobust beamforming scenario 

To verify the performance of the algebraic solution in the 
perfect-CSI case, we modeled a source node with M = 4 
antennas, 71 = 72 = 10 and pi = P2 = 10 W. The 
complex channel gain vectors were chosen randomly so that 
all components were complex Gaussian random variables with 
variance 1. We took (J2 = <Ji, and used a grid of value pairs 
(iTi, CTfl) covering the range 1 < ui, (7^ < 2. 5000 simulations 
were performed for each For each simulation, the 

semidefinite-programming solution was computed using the 
Matlab-based convex optimization system cvx [8], and cases 
for which the exact solution required power of more than 
25 W were discarded. For the remaining cases, the conjugate 
gradient algorithm was used to estimate the algebraic solution 
corresponding to the “+’ sign in (16). The starting point 
for the conjugate gradient was chosen as a: = y = 0, 
which corresponds to the “maximal-ratio receive, maximal- 
ratio transmit” (MRR-MRT) suboptimal solution in [3]. Iter¬ 
ation was terminated when the power reduction achieved by 
the latest iteration was less than 0.5 percent. 

Over the entire range of parameter values, the algebraic so¬ 
lution evaluated using conjugate gradient gave average power 
increases of less than 0.008 dB over the optimal semidefinite- 
programming solution. Convergence of the conjugate gradient 
solution required 2-5 iterations (on average) over the range of 
parameter values, as shown in Figure 1. These results confirm 
that t the solution of (16) with the ‘-f’ sign is optimal. 

The suboptimal MRR-MRT solution also performed very 
well, and produced power increases of only 0.05-0.35 dB as 
shown in Figure 2. Although the improvement of the conjugate 
gradient solution over the MRR-MRT solution is not great, 
it comes at very low cost; each conjugate gradient iteration 
requires only about 300 MAC operations, as compared to 
about 300 operations for the MRR-MRT solution itself. In 
contrast, the complexity of an exact solution via convex 
programming was computed in [1] as over 500,000 operations 
{0{k^ + n^k + n^k^) with k = 64, n = 8). 

VI. Simulation of robust beamforming scenario 

Reference [1] demonstrates that given a perfect-CSI solu¬ 
tion, a low-complexity suboptimal solution with very high per¬ 
formance can be found for the robust case with imperfect CSI. 


Fig. 1. Number of iterations until convergence for conjugate-gradient 
solution. 


Fig. 2. Transmit power increase over exact solution from MRR-MRT 
solution. 


In this section, we compare the performance of suboptimal 
robust solutions based on each of the three nonrobust solutions 
modeled in the previous section. 

Simulations were performed for a robust beamforming 
scenario with M = A, = IFF, cr| = IFF, ,Pi = lOFF, 
7i = 10, and = [0.01,0.15] with increments of 0.02 
{i = 1,2). The channel was generated as hi ~ CA/'(0,/); 
and the channel estimation error Ahi was generated as Ahi ^ 
CN'{0, (ef/16)/), which corresponds to a probability of 0.958 
that II Aliill < ei,i = 1, 2). 

In the simulations an outage was declared when the SINR 
at either source node fell below 7^. In Fig. (3), and (4), we 
















respectively plot the outage probability and 95*^ percentile 
of the empirical cumulative density function (cdf) [10] of the 
transmit power required to achieve the corresponding outage 
performance. The perturbed nonrobust solution described in 
[1] is used with each of the three nonrobust solutions: the 
semidefinite programming, algebraic, and MRR-MRT cases 
are indicated by “exact”, “conjugate gradient”, and “MRR- 
MRT” in the figures. As in the previous section, only one 
of the six candidate algebraic solutions was computed. There 
was virtually no difference between the performance of the 
exact and conjugate gradient solutions: this confirms the result 
of the previous section that one particular candidate of the 
six candidate solutions almost always gives the best overall 
solution. The simulation also shows that the conjugate gradient 
solution brings slight outage reductions (0-3%) and power 
reductions (1.5-4.5%) over the MRR-MRT approximation. 


Outage Probability- vs. Channel Estimation Error 



Fig. 3. TWR beamforming performance: outage vs. channel estimation error 
parameter. 


95th Percentile Relay Power vs. Channel Estimation Error 



Fig. 4. TWR beamforming performance: power vs. channel estimation error 
parameter. 

VII. Conclusions 

Simulation results show that the conjugate gradient imple¬ 
mentation effectively gives the exact solution to the beamform¬ 


ing problem (2). Furthermore, it may be used to construct 
a low-cost, high-performance solution to the corresponding 
robust problem. In the case of 2-way relays, the algebraic 
solution provides only limited performance improvements over 
the even lower-cost MRR-MRT solution. Nonetheless, our 
results demonstrate mathematical techniques for obtaining 
computationally inexpensive exact solutions to a nonconvex 
beamforming-optimization problem. Similar techniques may 
be applied to other situations to obtain lower-cost, higher¬ 
performing beamforming solutions. 
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